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\0 ■ ABSTRACT 

This paper generalises the hybrid power spectrum estimator developed in Efstathiou 
(2004a) to the estimation of polarization power spectra of the cosmic microwave back- 
ground radiation. The hybrid power spectrum estimator is unbiased and we show that 
it is close to optimal at all multipoles, provided the pixel noise satisfies certain rea- 
sonable constraints. Furthermore, the hybrid estimator is computationally fast and 
can easily be incorporated in a Monte-Carlo chain for Planck-sized data sets. Sim- 
ple formulae are given for the covariance matrices, including instrumental noise, and 
these are tested extensively against numerical simulations. We compare the behaviour 
^ ' of simple pseudo-CV estimates with maximum likelihood estimates at low multipoles. 

For realistic sky cuts, maximum likelihood estimates reduce very significantly the mix- 
ing of E and B modes. To achieve limits on the scalar-tensor ratio of r <C 0.1 from sky 
maps with realistic sky cuts, maximum likelihood methods, or pseudo-Cf estimators 
C^) ' based on unambiguous E and B modes, will be essential. 

\D '. 

f^) . Key words: Methods: data analysis, statistical; Cosmology: cosmic microwave back- 

ground, large-scale structure of Universe 
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■ In an earlier paper (Efstathiou 2004a, hereafter E04) a hybrid power spectrum estimator was developed for temperature 
| anisotropies of the cosmic microwave background (CMB) anisotropies. This estimator combined quadratic maximum likelihood 
' (QML) estimates at low multipoles with a set of pseudo-C^ (PCL) estimates at higher multipoles to produce a near optimal 

rS ' power spectrum estimate over the entire multipole range for realistic sky coverages, scanning patterns and instrument noise. 

■ We used analytic arguments and large numbers of numerical simulations to demonstrate that the method was near optimal 
and that the covariance matrix over the full range of multipoles could be estimated simply and accurately. 

A full maximum-likelihood estimator would require the inversion and multiplication of Nd x Nd matrices, where Nd is 
the size of the data vector. For WMAP (Bennett et al. 2003a) or Planck (Bersanelli et al. 1996; The Planck Consortia 2005) 
sized data sets, with Nd £ 10 6 — 10 7 , a brute force application of C(A|) methods is impossible computationally. This has 
been portrayed as a major computational challenge in CMB data analysis (e.g. Bond et al. 1999; Borrill 1999) and various 
approximate, but computationally demanding, techniques for solving this problem have been proposed (e.g. Oh, Spergel and 
Hinshaw 1999, Dore, Knox and Peel (2001), Pen 20003). One of the major motivations for E04 was to show that an 0(Aj) 
computation is unnecessary and that a near-optimal estimator (essentially indistinguishable from an exact 0(Aj) maximum- 
likelihood solution), with a calculable covariance matrix, could be estimated simply by combining PCL and QML estimates. 
Furthermore, the analysis of realistic experiments must deal with uncertainties such as correlated instrument noise, beam 
calibrations, foreground separation, point sources etc. It is simply not worth performing an exact C(A|) power spectrum 
analysis if the assumptions under which it is optimal are violated by 'real world' complexities, and if these complexities 
cannot be folded into the error estimates. With a fast hybrid estimator it is feasible to assess such 'real-world' complexities 
using a hybrid estimator within a Monte-Carlo chain and to quantify any corrrections to the covariance matrix, or likelihood 
function. The ability to analyse such 'real-world' effects is likely to be much more important than any hypothetical marginal 
improvement that might be gained by applying a full 0(A|) method. We refer the reader to E04 for a discussion of the 
rationale for a hybrid estimator. The arguments will not be repeated here. Instead, we will focus on generalising E04 to the 
estimation of power spectra of the CMB polarization anisotropies. 

As is well-known, the Thomson scattering of an anisotropic photon distribution leads to a small net linear polarization of 
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the CMB anisotropies (for an introductory review and references to earlier work see Hu and White 1997). This polarization 
signal can be decomposed into scalar B-modes and pseudo-scalar B-modes. The separation of a polarization pattern into E 
and B modes is of particular interest since scalar primordial perturbations generate only E mode while tensor perturbations 
generate E and B modes of roughly comparable amplitudes (Zaldarriaga and Seljak 1997; Kamionkowski, Kosowsky and 
Stebbins 1997). The detection of an intrinsic B-mode signal in the CMB would provide incontrovertible evidence that the 
Universe experienced an inflationary phase and would fix the energy scale of inflation (see e.g. Lyth 1984). 

An B-mode polarization signal was first discovered by DASI* (Kovac et al. 2002, Leitch et al. 2004). Exquisite measure- 
ments of the temperature- B-mode (TE) cross power spectrum have been reported by the WMApt team (Kogut et al. 2003). 
Measurements of the B-mode power spectrum have been reported by the OBI* experiment (Readhead et al. 2004) and by the 
2003 flight of Boomerang (Montroy et al. 2005). Primordial B-mode anisotropies have not yet been detected in the CMB. The 
detection of B-mode and possibly B-mode anistropies is one of the main science goals of the Planck mission and of ground 
based experiments such as Clover (Taylor et al. 2004). 

Associated with these experimental developments, there have been many investigations of techniques for analysing E and 
B modes from maps of the CMB sky. These analyses can be grouped, approximately, into the following catagories: 

(i) PCL estimators and correlation functions: Fast methods of estimating E and B-mode power spectra using correlation 
functions are described by Chon et al. (2004). Statistically equivalent PCL estimators are described by Kogut et al. (2003), 
Hansen and Gorski (2003), Challinor and Chon (2005, hereafter CC05) and Brown et al. (2005). In particular, CC05 present 
analytic approximations to the covariance matrices of PCL estimates for the case of noise-free data, while Brown et al. (2005) 
develop a Monte-Carlo method for calibrating covariance estimates from incomplete maps of the sky. PCL estimators can be 
evaluated using fast spherical transforms and hence scale as 0(N^ 2 ). 

(ii) Maximum Likelihood Methods: The generalization of the iterative maximum likelihood power spectrum estimation meth- 
ods of Bond et al. (1998) to the analysis of polarisation is straightforward and will not be discussed further here. Tegmark and 
de Oliveira-Costa (2001, hereafter TdOOl) define a quadratic estimator which is based on assumed forms for the temperature 
and polarization power spectra, generalising earlier work of Tegmark (1997). This method, which we will refer to as QML, is 
equivalent to a maximum likelihood solution if the guesses for the power spectra are close to their true values. As explained 
above, these methods involve matrix inversions and multiplications which scale as O(iVj). 

(iii) Harmonic E and B Mode Decomposition: The Stokes' parameters Q and U describing linear polarization define a rank 
two symmetric trace-free polarization tensor on sphere. Over the complete sky, the polarization tensor can be decomposed 
uniquely into E and B modes. However, since the E and B decomposition is non-local, it is non-unique in the presence of 
boundaries. In any realistic situation, a sky cut must be imposed to exclude contamination of the CMB signal by high levels 
of Galactic emission at low Galactic latitudes. Various authors have discussed ways of detecting pure E and B modes from 
Q and U maps on an incomplete sky (Lewis et al. 2002; Bunn et al. 2003; Bunn 2003; Lewis 2003). A key motivation for 
these analyses has been for diagnostic purposes (e.g. checking for a frequency dependent Galactic B-mode signal). However, 
as this article was nearing completion an interesting paper appeared by Smith (2005) describing how to construct pseudo-Cf 
estimators from unambiguous E and B modes on a cut sky. This type of analysis can effectively eliminate the severe E and 
B mode mixing that afflicts simple pseudo-CV estimators at low multipoles (see Section 3.4). (The method is not completely 
straightforward, nor unique, because it requires a 'pre-estimation' step to define unambiguous modes. We will throughout this 
paper use the abbreviation PCL to refer to the simple pseudo-C< estimators as defined in Section 2.) 

As we will show in this paper, these three analysis techniques are closely related. To illustrate how these methods are 
interelated it is useful to consider separately the cases of noise-free and noisy data: 

• Comparison of PCL and QML estimators for noise-free data: On a complete sky with noise-free data, PCL and QML 
estimates are identical. However, a sky cut will couple CMB modes over a range of multipoles AL, and this will lead to mixing 
of E and B pseudo- multipoles defined over the incomplete sky (see Section 2). Although PCL power spectrum estimators can 
be defined which give unbiased estimates of the true E and B mode power spectra in the mean, the mixing of E and B modes 
is reflected in large variances and couplings between the estimated power spectra at multipoles I £ AL (Section 2). A QML 
estimator applied to noise-free data unscrambles the E- and B-mode multipoles, in much the same way as the direct modal 
decompositions described by Lewis et al. (2002), Lewis (2003) and Smith (2005). The QML estimator returns almost optimal 
power spectrum estimates with smaller variances than a PCL estimator (Section 3). In particular, for realistic sky cuts an 
estimator that minimises E and B mode mixing, such as the QML estimator, is essential if one wants to probe low amplitude 
B-mode signals (tensor-scalar ratios r <JC 0.1). At high multipoles, t >• AL, E and B mode mixing becomes unimportant and, 
in the case of noise free data, QML and PCL estimators become statistically equivalent (Section 3). 

* Degree Angular Scale Interferometer 

t Wilkinson Microwave Anisotropy Probe 
t Cosmic Background Imager 
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• Comparison of PCL and QML estimators for noisy data: If the following conditions are satisfied: (a) the Q and U maps 
have identical noise properties; (b) the noise is uncorrected with variance per pixel of; (c) the estimates of E and B mode 
power spectra are for multipoles higher than some characterstic multipole Ln, where noise dominates; then one can show that 
an inverse-variance weighted PCL estimator is statistically equivalent to a QML estimator (Section 4.2). For the intermediate 
multipoles, Ln ijl AL, the optimal weighting for PCL estimates is intermediate between equal weight per pixel and inverse 
variance weighting. By combining a set of PCL estimates with different weights, it is possible to define a fast polarization 
estimator (in analogy with the temperature estimator discussed by E04) that is statistically indistinguishable from a maximum 
likelihood estimator over a wide range of multipoles. Dealing with strongly correlated noise is more problematic. Fortunately, 
for Planck-type scanning strategies, the noise pattern (in Q, U and T) should be accurately white at high multipoles (Efstathiou 
2005). A detailed correlated noise model should therefore only be required at low multipoles. It is straightforward to include a 
model for correlated noise in a QML estimate at low multipoles, provided the pixel noise covariance matrix can be estimated 
for low resolution maps (Section 4.2). 

The layout of this paper is similar to that for the temperature analysis presented in E04. PCL polarization estimates are 
discussed in Section 2 together with analytic approximations for covariances matrices in the absence of instrumental noise. 
The relationship of this work to the results presented in CC05 is discussed. Section 3 discusses QML polarization estimators 
for noise free data. A large set of numerical simulations is used to illustrate the effects of E and B mode mixing on an 
incomplete sky for both PCL and QML estimators. Section 4 discusses PCL and QML estimators including instrumental 
noise and Section 5 discusses a hybrid polarization power spectrum estimator. To keep the discussion simple, most of the 
analytic and numerical results presented in this paper refer to E and B mode power spectrum estimation. Our discussion 
of the temperature-i? mode cross power spectrum (denoted C x in this paper) is, intentionally, less complete since no new 
concepts are required for its analysis. Our conclusions are summarized in Section 6. 



2 ESTIMATION USING PSEUDO-CV 

We begin by relating the Q and U Stokes parameters (defined in direction n with respect to a spherical coordinate system 
eg, e^) and the spin two harmonics: 

Q(n) ± iU{h) = y^a± 2 <!m ±2Y em (n). (1) 

Note that the sign convention for polarization in this paper follows that of Zaldarriaga and Seljak (1997), which differs from 
the IAU convention (see Hamaker and Bregman 1996). The coefficients a±2i m are related to the E and B mode multipole 
coefficients af m and af m by 

a±2im — — (af m ± iafm)- (2) 

If we have an incomplete sky, we can compute 'pseudo-multipole' coefficients aj m , af m and af m by computing the following 
sums, 

oL = ^TiWiQiY; m (ei), (3a) 

i 

af m = -\ + iU)iWiQi iYi m + {Q- iU)iWitii -iYl m = ~ ^^{QiRtm + iU i RJ m )w i Q.i, (3b) 

i i 

af m = \ + iU)iWi^li 2Y e * m - (Q - iU)iWiQ.i -iY? m = % -'^{Q l RJ m + iUiR^WiQ.r, (3c) 

i i 

where in equations (3b) and (3c) 

Rtm ~ iYim + -2Ye m , Rim ~ lYim — -2Yi m , (4) 

f2j is the area of pixel i and Wi is an arbitrary weight function with spherical transform 

m m = y]wj{liYl m (6i). (5) 
From these pseudo-multipole coefficients, we can form the following PCL power spectrum estimates 

f^T _ 1 V~"^ ~T ~T* f,X _ 1 V~"^ ~T ~E* A£ _ 1 V""* ~E ~E* fiB _ 1 V""* ~B ~B* / c x 

~ (21 1) ' J a tm a lmi W — I 1) / ' a ^ mffl ^ m ' e ~ Jq£ 2. I a<?ma l"" 1 — I 1) / ' a lm a tm- (" ) 

m m m m 

The expectation values of these PCL estimates are related to the true values C T , C x , C E , C B by 
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(b) KpO mask 



(d) U 




Figure 1. An example of the input maps used in the simulations described in Section 2. The pixel size is 8 C = 1° and the CMB sky 
has been smoothed with a Gaussian beam of FWHM 9 B = 2° . The concordance ACDM model with the parameters given in the text has 
been assumed and a tensor component with r = 0.2 has been added to give an non-zero B-mode contribution. Figures (a), (c) and (d) 
show the T, Q, and U maps, while Figure (b) shows the WMAP KpO mask. 
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where the matrices M T , M x , etc are given by (Kogut et al. 2003) 
(2l 2 + 1) 



4tt 



•^(2£ 3 + l)W / , 



^3 



il i2 h 





(2£ 2 + 1) 
8tt 



^(24 + l)W, 3 (l+(-l) 2 



h (.2 £ 3 




li £2 4 
-2 2 



Mi 



Mi. 



M B ± f 2 — 



16tt 

(2£ 2 + 



+l)^ (2£ .,_ 1) ^. j(1+( _p'. .• ' * 

3 + 1)W, 3 (1 -(-!)-) 



16tt 







and We is the power spectrum of the weight function 



We = 



+ 1) 51 



(2£+l) 
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■■£i+£2+ 4, 



(8b) 
(8c) 



(9) 



(cf equation (5)). In the absence of parity violating physics in the early universe, the primordial C TB and C EB spectra 
should be identically zero. Even if there is no compelling motivation from fundamental physics, there may be other reasons for 
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Figure 2. PCL polarization power spectra. The blue points show the convolved power spectra, Ce, averaged over 10 5 simulations with 
the same parameters as those in Figure 1 and including the KpO mask. The green points show the deconvolved spectra Ce . The solid lines 
show the forms of Cg and Ce computed from the theoretical input spectra. The dashed lines in Figures (a)— (c) show the contribution to 
Ce from the tensor component. 



wanting to estimate these spectra, e.g. for consistency checks and to test for systematic errors. However, to keep the analysis 
as simple as possible, they will be not be included in this paper. 

If the sky cut is small, then the matrix M appearing in equation (7) will be non-singular and hence invertible. If this is 
the case, then unbiased estimates of the true power spectra can be formed from the PCL estimates by evaluating § 

Ct = M-}C e ,. (10) 

As an example, Figure 1 shows noise-free T, Q and U maps from a single Gaussian realisation of the concordance ACDM 
model favoured by WMAP (Spergel et al. 2003). (The exact parameters adopted in this paper are those of the ACDM model 
defined in Section 2 of Efstathiou (2003)). The maps have been created using software written by the author that uses an 
'igloo' pixelisation scheme. The maps shown in Figure 1 have a pixel size of 1° and a symmetrical Gaussian beam smoothing 
of FWMH of 9s — 2° . In these simulations, we have assumed a tensor-scalar ratio of r = 0.2, where r is defined in terms of 
the relative amplitudes of the ensemble averages of the temperature power spectra at £ = 10, 

C T tensor 
f~l 1 eaalar 

As in E04, unless otherwise stated, beam functions will not be written explicitly in equations and so Ce will sometimes mean 
Cebj, where be is the Gaussian beam function 

be = exp + 1)(0.425<9 S ) 2 ) . (12) 

Figure 1 also shows one of the Kpn family of WMAP masks (see Bennett et al., 2003b, for a discussion of the WMAP 
masks). The KpO mask shown in Figure 1 removes about 21 per cent of the sky at low Galactic latitudes and is a relatively 
conservative Galactic mask. (For reference, the WMAP TE analysis reported by Kogut et al. (2003) used the less conservative 
Kp2 mask which removes around 13 per cent of the sky). 

Figure 2 shows the averages of PCL estimates Ce and Ce computed from 10 5 noise-free simulations with the same 
cosmological parameters used to construct Figure 1 and with the KpO sky mask applied. For the KpO mask, the matrix M in 

§ The notation here follows the notation introduced in E04. The expectation values of the hatted PCL estimates (fig) (if they are well 
defined) are equal to the true power spectra. The expectation values of the tilde estimates (fit) are given by convolutions of the true 
power spectra (equation 7). 
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KpO unapodised 



KpO apodised 




I I 

Figure 3. The diagonal components of the PCL power spectrum covariance matrix estimated from numerical simulations compared to 
the theoretical dispersions given in equations (15b) - (15d). The diagrams to the left show results for the unapodised KpO mask and 
those to the right show results for the apodised KpO mask illustrated in Figure 4. 



equation (10) is non-singular and can be inverted. Figure 2 shows that the PCL estimates Ci provide unbiased estimates of 
the true power spectra. 

For the temperature polarization power spectrum, accurate expressions for the covariance matrices of PCL power spec- 
trum at high multipoles (I ^> AL) can be derived quite easily (see E04). The analogous problem for polarization esti- 
mates is much more difficult (see CC05) because the covariances depend on the matrix products ±I(t m )(LM) ±I(LM)(t'm r ), 
±I(im)(LM) i^(iM)(fm'), where 

±I(tm)(t' m ') = j ^nw(n)( 2>7 m (n) 2 lW(n) ± _ 2 F/ m (n) - 2 Ye'm'(n)). (13) 

CC05 show that these products can be expressed as integrals ±Ium)(£'m, 1 ) as in equation (13), but with the window finction 
w replaced by w and integrals involving the gradients of w. For PCL estimates from noise free data, the covariance matrices 
can therefore be approximated by terms which depend on the coupling matrices in equations (8a - 8d), but with We replaced 
by the power spectrum of the square of the window function, 

we m = ^2 w1VliYe* m (ft) , (14) 

together with complicated terms that depend on gradients of the window functions Wi^ Ignoring the gradient terms, the 
covariances can be approximated by 

<ACfAC?>«M^MS„ (15a) 



^ There arc no such gradient terms for the temperature covariance matrix (ACjACj,), see E04 for a detailed analysis. 
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KpO mask apodised 




Figure 4. Apodised KpO mask at a pixel resolution of 6 C = 1° created using the iterative algorithm described in the text. 



(ACf ACf,) - Mf v + (|rf%^. (15b) 

(AC?ACS)h0$LmZ?, (15d) 

(ACf A(5£> - Cf ' )X '^}^' m (15e) 
The covariance matrix of the PCL estimates Ci is given by 

(AC e AC e ,) =M- 1 (AC i ACt,)(M~ 1 ) T , (16) 

where M is the matrix appearing in equation (7) (i.e. computed from We, rather than from the power spectrum of vjf). 
Neglect of the gradient terms will be referred to somewhat loosely as the 'scalar approximation', since the expressions (15a- 
15e) depend on only the spin-0 multipoles of the square of the window function (14). It is interesting to ask under what 
circumstances the scalar approximation provides accurate estimates of the covariances. For relatively small sky cuts, such as 
the KpO mask, these expressions provide quite accurate estimates for the X and E power spectra at multipoles £ AL, where 
AL is the mode-coupling scale introduced by the sky cut. The left-hand panels in Figure 3 show the diagonal components of 
the covariance matrices for the C? , C? and C? power spectra estimated from equations (15b) - (15d) compared to estimates 
from the 10 s simulations used to generate Figure 2. For the KpO mask, the characterstic coupling scale is AL ~ 20, and one 
can see from Figure 3 that the simple analytic expressions provide quite accurate estimates of the errors of the X and E 
power spectra at multipoles £ ^ 50, but fail by a factor of ~ 2 or more for the B component. The failure of simple analytic 
expression for the B-component is related to mixing of the E and B modes. For the E power spectrum, mixing of E and B 
modes is always unimportant at high multipoles, £ 3> AL (since the _E-mode amplitude is fixed by the dominant scalar mode). 
We therefore expect the scalar approximation to provide accurate estimates of the covariances for the X and E power spectra 
at high multipoles. However, if the intrinsic amplitude of the B mode is low, mixing of E and B modes can dominate the 
estimates of C? and the covariances {AC? AC*) and {AC? AC?) at high multipoles. For our chosen tensor-scalar amplitude, 
r — 0.2, mode-mixing dominates the C? amplitude at high multipoles, which is why the scalar approximation fails so badly. 

Intuitively, one would expect that it would be possible to significantly reduce the effects of E and B mode mixing at 
high multipoles by apodising the sky cut. Apodisation down-weights data close to the edge where the non-local nature of the 
E and B mode decomposition causes ambiguity between modes (c/CC05). Thus we would expect the scalar approximation 
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Figure 5. Covariance matrices for the B-mode polarization spectra, illustrating the effects of apodization on the PCL estimator. Figures 
5(a) and 5(c) show the covariance matrices computed from 10 5 simulations with the same parameters as those in Figure 1. The simulations 
used in Figure 5(a) use the unapodised KpO mask while those used in Figure 5(c) use the apodised KpO mask shown in Figure 4. Figures 
5(b) and 5(d) show the corresponding analytic covariance matrices using the approximation of equation (15d) 



to be much more accurate if we apodise the KpO mask. Of course, in this context apodisation will lead to some information 
loss, but in practice the purpose of applying a mask is to reduce the effects of systematic errors caused by inaccuracies in 
modelling Galactic emission. Hence, the data close to the edge of a mask will usually be less reliable than data far removed 
from the mask. The slight loss of information caused by apodising is therefore likely to be more than compensated by the 
reduction in E and B mode mixing at high multipoles. 

Figure 4 shows an example of an apodisation algorithm applied to the KpO mask. The apodisation algorithm works 
iteratively as follows. First, the KpO mask in this pixelisation is specified by Wi = or 1 according to whether pixel i lies 
within or outside the mask. Then: 

(i) compute the spherical transform of the set of weights, wp m , (equation 5) and multiply these coefficients by a smoothing 
factor of bt (equation 12) with Gaussian FWHM of S = 10°; 

(ii) construct a smoothed weight function Wi by performing an inverse spherical transform; 

(iii) set the weights Wi that lie within the original KpO mask equal to zero; 
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(iv) go back to step (i) and transform the smoothed weights to make a new set of coefficients wi m and continue steps (i) - 
(iv) for a specified set of iterations. 

The example shown in Figure 4 shows the apodised mask after 4 iterations. The effective smoothing length therefore has 
a Gaussian FWHM of 9 S = 20° and one can see that the 'islands' interior to the masked region evident in Figure lb have been 
smoothed to low amplitudes. The right-hand panels of Figure 3 show the diagonal components of the covariance matrices for 
the X, E and B mode PCL power spectra of the apodised maps compared to the predictions of the scalar approximation. 
The X and E dispersions now agree almost perfectly with the scalar approximation for I £ 30. For the B mode, the variance 
caused by E and B mode mixing at I £ AL is much less than for the unapodised case shown in Figure 3. For the diagonal 
components of the B mode power spectrum, the scalar approximation is accurate to a few percent or so. Figure 5 shows 
the structure of the B mode covariance matrices for the unapodised and the apodised KpO masks. In the apodised case, the 
covariance matrix is diagonally dominant at £ >, 20 and the diagonal components are accurately described by equation (15d). 

The scalar approximation is therefore a very good approximation for all three power spectra X, E and B. Of course, 
even with apodisation, the scalar approximation for the B mode spectrum will fail at some critical value of the tensor-scalar 
ratio r, hence one cannot simply rescale equation (15d) by the appropriate amplitude of Cf for arbitrarily low values of r. 
For low values of r, one can explicitly include the gradient terms which depend on the amplitude of CfCf, (equation (80) 
of CC05) to get an accurate model of the covariance matrix at high multipoles. However, for most forseable experiments, 
including Planck, instrument noise is likely to dominate the B-mode covariance matrix at high multipoles. In Section 4 we 
will show that in this situation the scalar approximation is extremely accurate independent of the intrinsic amplitude of Cf . 



3 QUADRATIC MAXIMUM LIKELIHOOD 
3.1 Preliminaries 

A QML estimator for temperature CMB power spectra is discussed by Tegmark (f997). The generalization of this estimator 
to the estimation of CMB polarization is straightforward and is discussed in detail by TdOOf . The estimator will be reviewed 
here briefly. We work in pixel space and define an input data vector x consisting of the temperature differences and Stokes 
parameters Q and U (defined with respect to a fixed coordinate system, as in Figure f ) specified at each pixel, 

x=(AT,Q,U). (17) 

The optimal QML power spectrum estimate is (TdOOl) 

y r t = XiXjETf, r = (T,X,E,B), (18) 
where the matrices E rt are given by 

and C is the covariance matrix of the data vector x, 

qTT qTQ qTV 

Ca = {x iXj ) = ( C QT C QQ C QU | . (20) 

qUT qUQ qUU 

The matrices E rl defined in equation (19) give formally minimum variance power spectrum estimates if the covariance matrix 
dj is set equal to the true covariance matrix. However, as noted by TdOOl, equation (18) mixes AT components of the data 
vector x with Q and U components in providing estimates of the E and B mode power spectra. This is undesirable because 
for realistic noisy data, systematic errors in the AT measurements could contaminate estimates of the much lower amplitude 
E and B mode power spectra. It is therefore safer to separate the AT measurements from the Q and U measurements by 
'reshaping' the matrix (20): 

C TT 

<5 I3 = | C QQ C QU I (21) 

C UQ C uu 

and using the matrices 

2 dC\ y ' 

in place of the matrices E rl in equation (18). The estimates y r e will give unbiased estimates of the true power spectra Cf , 
s = (T,X,E,B), 
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(y r i) = F e s [,C(,, 

where 



dc dc 



dCl 



dc: 



(23a) 



(23b) 



By reshaping the covariance matrix, the estimate (18) will no longer be minimum variance. However, it is straightforward to 
show that for a KpO-type mask, the increase in variance caused by using (21) instead of (20) is almost imperceptible. 
Following E04 it is useful to define re-scaled power spectra 



(24) 



which have similar shapes to the true power spectra C\. Furthermore, if the matrix F is invertible, one can define unbiased 
estimates of the true power spectra via 



O'i = b y. 

These estimators can be considered the QML analogues to the PCL estimators defined by equations (6) and (10). 
The covariance matrices of the QML estimates is given by 



CE rt CE st ' 



(25) 



(26) 



where F rs is the Fisher matrix. From equation (25), the covariance matrix of the deconvolved QML estimates is given by 

(27) 



(ACVACV) = F~ 1 FF~ 1 



Notice that if we had used the true pixel covariance matrix (20), rather than the 'reshaped' form (21), the covariance matrix 
(26) would take the more familiar form 



F% = ^Tr 



dc c _! dc c _! 



dCl 



Notice also from equations (18) and (22) that the estimator y\ can be written as 
_ (21 + 1) 



20. 2 



Cf 



where the C zr are the PCL power spectra as denned in equation (6) but computed from the maps 

Zj = Cji Xi. 



(28) 



(29) 



(30) 



(fi in equation (29) is the solid angle of a single map pixel, all assumed to be of identical area) Thus, if one can invert 
the matrix C 1 ^ 1 , the QML estimator can be computed rapidly by applying fast spherical transforms to the maps defined 
by equation (30). This is particularly useful if one wants to compute Monte-Carlo simulations of the QML estimator to 
characterise the covariance matrix. For applications to maps with large numbers of pixels, it may well be more efficient to use 
Monte-Carlo simulations to estimate errors, rather than brute force evaluation of (26), particularly if the covariance matrix 
(27) is accurately band diagonal (which is true for KpO-like sky masks, see Figure 7 below) . 



3.2 Covariance Matrix at High Multipoles 

In this Section, we will assume that Cf = and derive an approximation to the covariance matrix of Cf for high multipoles 
and a small sky-cut. The analysis presented here is similar to the analysis of the QML temperature power spectrum covariance 
matrix given in Section 3.4 of E04. 

If the B-mo&e is absent, all of the information on the polarization anisotropies is contained in the spin 2 field 2-Pi = 
(Qi + iUi). The spin —2 field (Q; — iUi) contains no additional information. We can therefore construct a quadratic estimator 

yf= iPiWiiPjWjEa, (31) 

where Eij is as defined in equation (19) with, 

Cij = ( 2 PiWi 2-P/wj) = ^2 CfwiWj 2 Y em (i) 2Y e * m (j). (32) 

where Wi is an arbitrary pixel weight function. If the sky cut is small, it will be a good approximation to set the inverse of C 
to 
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^-E^^-^^O') (33) 
using the completeness and othogonality relations for the tensorial harmonics. Hence, 

^ffe * E aY ^ w 2YL{k) (34) 

and the Fisher matrix (23b) is 

Fu> = \c^^c^ « 2 2 an, 2 y £m (i) 2 Y; m (k) 2 Y e , m , (k) 2 r;, m , (i) 

*~ ^ ^ ^ ik mm' 

- W + V-MS? (35) 



ns~iE s~iE 
ZO <? ° ell' 



where M EE is the coupling matrix of equation (8c). Equation (34) agrees with equation (15c) and shows that for small 
sky-cuts and noise-free data, the QML estimator at high multipoles (I 3> AL) is statistically equivalent to the PCL estimator 
with equal weight per pixel. 

The derivation given here is only valid in so far as the B-mode can be ignored and the sky cut is small enough (and the 
multipoles high enough) that the approximation (33) is accurate. It is not straightforward to generalise the analysis presented 
here to admixtures of E and B modes. Given the results of Section 2, we would expect that at high multipoles and low 
B-mode amplitudes, the QML estimator for the B-mode power spectrum will be statistically equivalent to a PCL estimator 
with suitably apodised sky masks. In other words, the QML estimator will self-consistently find the appropriate apodisation 
(via equation 30) for any given geometry to minimise E and B mode leakage at each multipole. 



3.3 Covariance Matrix at Low Multipoles 

Generalising the discussion in Section 3.3 of E04 to polarization, we can consider the vector Xi to consist of the harmonic 
coefficients af m and af m measured from noise free data on the cut sky. These are related to the true harmonic coefficients af m 
and af m by a coupling matrix 



O-fm \ ( ~ +I{tm){l'm') ~ * -I(tm)(t'm') \ ( O.'j 

-i -I(lm)(i'm') +I(lm)(l'm') 



(36) 



where the integrals ±I(t m )(i'm') are defined in equation (13). Evidently, if the coupling matrix in equation (36) can be 
inverted, then the true power spectra can be reconstructed. However, if the power spectra extend to high values of t and 
m, then the coupling matrix will be singular for a finite sky cut, since it is impossible to reconstruct modes that lie entirely 
within the cut. Furthermore, the presence of boundaries leads to 'ambiguous' modes, i.e. modes which cannot be assigned as 
pure E and B modes. As discussed in the Introduction, various authors have investigated harmonic mode reconstruction from 
data on an incomplete sky (Lewis et al. 2002; Bunn et al. 2003; Lewis 2003). The analysis by Lewis (2003) is particularly 
relevant to the discussion because he shows how to construct projections of pure E and B modes for noise-free data with a 
Kp2 sky mask. For low values of I, {I H l/8 cu t, where 8 cut is the characteristic width of the sky cut), the true E and B modes 
can be reconstructed accurately, (i.e. for suitably band- limited data the coupling matrix (36) is invertible. As in the case 
of temperature anisotropies discussed by E04, for noise-free data on a cut-sky the QML estimator will effectively invert the 
coupling matrix at low multipoles, reconstructing almost the true values of the E and B mode power spectra for the full sky. 
The covariance matrix of the QML estimators at low multipoles should therefore be given by cosmic variance, i.e. diagonal 
with components 

At high multipoles, the variance of the QML estimates will increase above cosmic variance (as given by equation (35) for the 
E-mode) to 

where / s k y is the fraction of the sky surveyed, reflecting the loss of information on modes within the sky cut. (See CC05 for 
the analogue of equation (38) for {(ACf) 2 } at high multipoles in the limit Cf = 0.) The expectations of equations (37) and 
(38) are borne out by the tests described in Sections 3.4 and 5. 
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Figure 6. Comparision of PCL and QML polarization power spectrum estimates using the simulations described in Section 3.4. These 
simulations are noise-free, smoothed with a Gaussian of FWHM 9 S = 7° and have had the KpO mask applied. In the case of the PCL 
estimates, an apodised KpO mask was applied using the algorithm described in Section 2 (see Figure 4). The (blue) stars in each panel 
show the convolved estimates Ci and the (red) circles show the deconvolved estimates C'{. The (green) lines in Figures (a) and (c) show 
the Cf spectrum of the ACDM model. The input Cf spectrum was assumed to be zero. Error bars show la errors. Note the large 
difference in scales on the diagrams showing the PCL and QML B-mode power spectrum estimates. 



3.4 Comparison of PCL and QML methods for the KpO mask 

Using the 'reshaped' covariance matrix of equation (21), the problem of estimating E and B mode power spectra becomes 
decoupled from the problem of estimating T and X power spectra. We will therefore confine the discussion in this Section to 
the estimation of E and B mode power spectra using the QML methodology described in Section 3.1. To compare the QML 
and PCL methods we have performed a set of 10 5 noise-free simulations of the concordance ACDM model with 5° pixels and 
a smoothing of 9 S — 7° . The KpO mask was applied, leaving 1216 active pixels out of a total of 1632 pixels over the full sky. 
For these simulations, the amplitude of the input B mode power spectrum was set to zero. For the PCL estimates, the KpO 
mask was apodised as described in Section 2 (see Figure 4). 

Results from these simulations are shown in Figures 6 and 7. Figure 6 shows that the deconvolved PCL and QML 
estimates provide unbiased estimates of the input power spectra. The errors of the E-mode spectra are dominated by cosmic 
variance, though one can see that the KpO sky cut causes a noticeable increase in the PCL errors. This is equivalent to the 
discussion of 'estimator induced' variance for the temperature anisotropies at low multipoles discussed in Section 3.3 of E04 
(see also Efstathiou 2004b). At low multipoles (£ 10) the QML estimator returns very nearly the exact value of Cf for each 
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(a) simulations (PCL) (b) analytic (PCL) 






Figure 7. Comparison of the covariance matrices for PCL and QML estimates of the E and B mode polarization spectra. The upper 
diagrams show the covariance matrices of the 'deconvolved' PCL estimates Cf and CP (equation 10) derived from numerical simulations 
(see Figure 7) compared with the analytic predictions of equations (15c)- (15e). The lower diagrams show equivalent plots for the QML 
estimates Ci (equation 25) compared with the analytic estimates derived from equation (27). (Note that the EB quadrants of these 
diagrams have been rescaled by large factors to render them visible.) 



realisation and so the covariance matrix is accurately diagonal (Figures 7c and 7d) and limited by cosmic variance over the 
full sky (equation 37). In contrast, the PCL estimates of Cf at low multipoles are strongly correlated (Figure 7a) as a result 
of the sky cut and have a large estimator-induced variance which amplifies the errors above those of the QML estimator. 

The contast between the PCL and QML B-mode estimates is even more extreme. The input B-mode spectrum in these 
simulations is set to zero, thus any measured B-mode is arising from E mode leakage caused by the sky cut. From Figures 
6b one can see that for the PCL estimator, leakage leads to B-mode amplitudes of ~ 0.01-0.02 (fiK) 2 , i.e. about 10% of the 
amplitude of the B-mode spectrum, simply as a result of the sky cut. The estimator induced variance is therefore huge and 
so even if one had noise-free data, it would not be possible to probe tensor-scalar ratios of r <C 0.1 using a KpO sky cut and 
PCL power spectrum estimates (see also CC05 for a similar conclusion). As expected from the discussion in Section 3.3, the 
QML estimator recovers almost the exact input power spectrum at low multipoles even on the cut sky. Thus in Figure 6d 
the QML B-mode estimates have an amplitude of typically 10 -4 ()J-K) 2 and the covariance matrix for these estimates is 
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accurately diagonal (Figure 7d). The residual 'estimator induced' variance is a consequence of the coupling between low I 
modes and high t modes which limits the accuracy of an inversion of equation (36). Thus, for a KpO-type sky cut the variance 
of the QML estimator leads to a limit on the detectability of a tensor-scalar ratio of r ~ 1CP 5 in the absence of noise or 
residual foreground emission (see also Amarie, Hirata and Seljak 2005). It would be possible to probe lower values of r by 
isolating pure B-modes on the cut sky, as described by Lewis (2003), though methods for probing r ,$ 10~ 5 are probably only 
of academic interest at this stage, given our lack of knowledge of polarised foregrounds. 



4 INCLUDING INSTRUMENTAL NOISE 
4.1 PCL estimator with noise 

In this Section we consider the effects of simple models of instrument noise. The data vector (17) is therefore split into a 
signal contribution x" and a noise contribution x n , which are assumed to be uncorrelated with covariance matrices S and N: 

Xi = xt + x?, Si, = (x 3 t x°), Nij = (39) 

We assume further that the noise on each of the T, Q and U maps is uncorrelated. With these assumptions, the expectation 
values of the PCL power spectrum estimates is modified from (7) to 

(Cf) = M T C T + Nj, (40a) 

(Cf) = M x C x , (40b) 

{Cf > = M EE C E + M EB C B + Nf, (40c) 

{Cf ) = M BB C B + M EB C E + N B , (40d) 
where the noise power spectra Nj , N E and N B are given by 

Nj = (2(^+1) Y trn{i)Ye.rn{j)WiWjQ.iSl j , (41a) 
ijm 

NF = p^ijjE^^W^^ + N^R^m^U))^^, (41b) 

ijm 

NfB = (aTi)jE^ E «™W^(^ + ^Rt*(t)RL(j))™^^. (4ic) 

ijm 

If the noise covariance matrices are diagonal, 

W« = (ffD 2 *«, N^ = {aff8 %J , A^ = (a?) 2 5 l3 , (42) 
then the expressions (41a)-(41c) simplify to 

W/ = ^5>nV«?, (42a) 

i 

= = ^ J>«) 2 + (aYf)wfrl (42b) 

i 

where the last equation follows from the addition relation for tensor ial harmonics, 

s.YeUOiAi) s 2 Y im (6 2 A2) = (^-) 1/2 s 2 Y e - ai ([3,a)e-^\ (43) 

m 

and the angles (#i,0i) and (02,</>2) are related to the Euler angles (a,/3,j) as described in Section 1.4.7 of Varshalovich, 
Moskalev and Khersonskii (1988). 
Defining the coupling matrices, 

HH^,*) = £^*,,( f ' \ *)'. (44a, 
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ZMlul*,W) = J2^^Wz 3 (l + (-l) L f( % "l "l \ , (44c) 

a Bfl( / ll / a ,wO = E^^wi,(i-(-i) i ) a (^ l l *Y, (44d) 

for arbitrary W, the analogues of equations (15a) - (15e) including uncorrelated instrumental noise in the case {of) 2 = {oY) 2 
are: 

{AC? AC?,) « 2CJ , QTHt(£,/,# 2 ) + 2Ht(^/,# tt )+4(Q t C , ^) 1/2 H t (£,/,# 2T ), (45a) 

{AC? AC?;) « {C?C?,C?C?) 1/2 E X {£,£', W 2 ) + C?C?E T {£,£', W 2 ) + E x (£,£', W TQ ) + 

{C?C?,) 1/2 E X {£,£',W 2Q ) + {C?C?,) 1/2 E X (£,£',W 2T ), (45b) 

{AC? AC?) w 2CfCj?H i5i5 (£,^W' 2 ) + 2H BB (^ (45c) 

{AC? AC?,) w 2CfC^SBB(^^,W r2 ) + 2SBB(€,^,W 00 ) + 4(CfC^) 1/2 SBB(<,^,W r2 °), (45d) 

{AC? AC?,) « [{C?C?) 1/2 + {C?C?) 1/2 ] 2 E eb {£J\W 2 ) + 2Eeb{£J\W qq ) + 4{C?C?,) 1/2 Eeb{£J\W 2Q ) (45e) 
where the various window functions appearing in equations (45a)-(45e) are defined by 

T = T^yy J! ^( W L»;™), = J4> Tl7 = — - £ &(»L^), (46a) 



2 



^ = (2£TT) E ^M)' ^ - ^ - (2£Tl) E ( 46b ) 

m m 

= W QU = W uu = ]T K Q J 2 , (46c) 

m 

and 

"4m = ^ w 2 f2;Yft„(j), Wfai = y~](crT ) 2 w 2 Q 2 Ye m {i), wf m = w u im = ~^^{cyf) 2 w 2 ^l 2 Yi m {i) . (47) 

i i i 

To check these expressions, we repeated the simulations described in Section 2 {i.e. 9 C = 1°, 6 S = 2°, r = 0.2) but 
including uncorrelated noise with cr q = a u = cr T /\/2 = 2.725^K. With these noise levels, the B-mode power spectrum is 
noise dominated for I £ 20. PCL power spectrum estimates, using the apodised KpO mask of Figure 4, were computed for 10 5 
simulations. 

Figure 8 compares the diagonal components of the covariance matrices for the X, E and B power spectra compared to the 
analytic expressions of equations (45b) - (45d). Apart from the addition of instrumental noise, these simulations are identical 
to those described in Section 2 and so Figure 8 can be compared directly with Figure 3. Figure 9 compares the full forms of 
the covariance matrices computed from the simulations to the analytic predictions of equations (45b) - (45e). Figures 8 and 9 
show that that the scalar approximation is extremely accurate for all three spectra X, E and B for £ £ 20. In the discussion 
of Figure 3, we pointed out that in the noise free case, the scalar approximation would break down for low tensor amplitudes 
(when the PCL B-mode amplitude is dominated by E-B mode mixing rather than true B modes). However, if the B-mode 
power spectrum is noise dominated at high multipoles (as it is in the simulations described in this Section), the covariance 
matrix is fixed by the pure noise term in equation (45d) which is independent of the amplitude of the intrinsic B mode. Thus 
in the noise dominated case, the scalar approximation will provide an accurate estimate of the PCL covariance matrices at 
high multipoles independent of the amplitude of the tensor component. 



4.2 QML estimator with noise 

Evidently, if the data vector x% includes instrumental noise, the generalization of equation (18) is, 

y\ = XiXjEll - N i:j E r 4, (48) 
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KpO apodised+noise 




I 

Figure 8. The diagonal components of the PCL power spectrum covariance matrix estimated from numerical simulations including 
instrumental noise (Section 4.1) compared to the theoretical dispersions given by equations (45b) - (45d). Apart from the addition of 
noise, the parameters of the simulations were identical to those used to generate Figure 3. 



and the covariance matrices C in equations (19) and (20) are replaced by dj = Sij + Nij, where S and N are the signal and 
noise covariance matrices defined in equation (39). 

If the noise is diagonal, N[f = (ol) 2 8ij8 rs , then in the noise dominated limit C^ 1 ~ 1/crfSij and one can see immediately 
from equation (30) that the QML estimator is mathematically equivalent to a PCL estimator applied to a map with inverse 
variance weighting, w[ = l/(a r ) 2 . As in the discussion of temperature power spectrum estimation (Hinshaw et al. 2003, E04), 
we find that the optimal weighting of a PCL estimator is: 

(i) equal weight per pixel in the signal dominated limit; 

(ii) inverse variance weighting in the noise-dominated limit, if the noise covariance matrix is diagonal. 

As described in E04, one of the key motivations for a hybrid estimator is to construct an estimator that combines 
PCL power spectra computed with various weighting schemes to give a near optimal power spectrum estimate between the 
two extremes (i) and (ii). This intermediate regime is difficult to analyse analytically. The hybrid estimator is discussed 
in the next Section, which parallels closely the discussion for temperature anisotropics given in E04. There is, however, an 
important difference between the discussion of temperature and polarisation anisotropics. Consider, for example, the optimal 
inverse-variance weighted i?-mode pseudo-multipole, 

5 - = ~\ E ( -$y?& + ' (jfj^) * (49) 

Because of the spin ±2 nature of the polarisation anisotropies, it is not possible to derive simple expressions analogous to 
equation (7) relating PCL power spectra constructed from these multipoles to the true power spectra unless {of ) 2 = (oY) 2 ■ If 
(°i ) 2 7^ (°i0 2 j ono would need to evaluate (time-consuming) products of t 2 x I 2 matrices to recover the true power spectra. 
For simplicity, we therefore assume (of) 2 = (erf 7 ) 2 in the rest of this paper, recognising that if this is not true, then it is not 
possible to define an easily computable optimal PCL estimator in the noise-dominated limit. 
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(a) simulations X— X 



(b) analytic X-X 




(c) simulations E— E 



(d) analytic E-E 





Figure 9. The PCL X — X, E — E, B — B and E — B covariance matrices for the simulations including instrumental noise (Section 4.1) 
compared to the analytic approximations of equations (45b)-(45e). 



5 HYBRID POLARIZATION POWER SPECTRUM ESTIMATOR 



Following E04, if we estimate a number of power spectra C\ a derived from the same data, where r denotes the mode 
(r = T, X, E, B) and a denotes the estimator (e.g. PCL estimates with different weight functions), then they can be combined 
into a single data vector C^i from which we can define a x 2 , 



(50) 



where Cf h is the hybrid estimator for mode r, T^t.pu is the inverse of the covariance matrix (AC^ AC%i ). Minimising 
equation (50) gives the following linear equations, which we can solve to form the hybrid estimate C( h 



•r al-LPl^t., - 7 

with covariance matrix 



sa/3i 2 



(51) 



(AC^AC^) 



7 , 3 r a i 1 pe 2 



(52) 



Evidently, since the hybrid estimator involves linear combinations of the power spectrum estimates (7J, it will provide an 
unbiased estimate of the true power spectra, provided that each of the CI is unbiased. Nevertheless, as with the discussion 
of the optimal QML estimator in Section 3.1, the solution of (51) mixes power spectrum estimates for different modes, which 
may be undesirable in practice. For example, for PCL estimates, the covariance matrix {ACf AC*} is zero by symmetry, thus 
the PCL estimates Cj and C* are strictly independent. However, the covariance matrices {ACf ACfj) and {ACj ACf,) are 
non-zero and are given by 



(ACjACp) 



(2£+l)(2£' + l) 



E 



E 

£ i mi 



K(l m )(l imi ) +I(t' m ')(l 1 m 1 ) 



(53a) 
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I I 

Figure 10. The points show the hybrid estimator applied to PCL estimate of the E and B mode power spectra as described in the text. 
The dashed lines show the white noise level. The lines show the power spectrum of the fiducial ACDM model and the error ranges given 
by the diagonal components of the covariance matrix (56). As the B mode power spectrum is noise dominated for most of the multipole 
range, the points with error bars in Figure 10b show the averages and dispersions of the power spectrum estimates in bands of width 
A£ = 50. The theoretical error ranges in this figure have been scaled appropriately for bands of width A£ = 50. 



(ACjACf,) 



(2£+l)(2l' + l) 



' + D ^ 



K(lm){tim l ) -1(1 



(e'm')(e 1 m 1 ) 



(53b) 



where ±I(e m )(e'm') are the coupling integrals denned in equation (13) (with w set to the weight function applied to the Q 
and U maps) and ^(( m )(i' m 'j is the coupling integral 

(53c) 



#"(£m)(£'m') = / dnw(n)Y e * m (h)Y e/rn ,(n), 

with w set to the weight function applied to the temperature map. The calculation of the full solution of equation (51) is 
therefore complicated, and unless the errors are very well understood, there is a danger of introducing biases by mixing modes 
with very different amplitudes. Rather than use the strict minimum variance solution (51), it is preferable to ignore cross 
covariances and to compute for each mode separately the solution of 



a/3l 2 



(54) 




1 e"si'C 1 i 



(55) 



The estimates (55) will be technically sub-optimal, but for practical purposes they will be indistinguishably close to optimal. 
The covariances and cross-covariances of the hybrid estimates are straightforward to compute and are given by 



B r i t a l (ACatACpti ) B 



(56) 



To illustrate the method, we generated a simulated sky of the concordance ACDM cosmology with r = 0.2, at a pixel 
size of 8 C = 0.1° and with a Gaussian smoothing of 9 S = 0.25°. We simulated the hit-count distribution for a Planck- type 
scanning strategy as described in E04, i.e. a single detector pointing at 85° to the spin axis of the spacecraft, and with the 
spin axis precessing slowly about the ecliptic plane according to 5°sin(20 e ), where 4> e is the ecliptic longitude. With this type 
of scanning strategy, regions with high hit counts are concentrated at the ecliptic poles as shown in Figure 11 of E04. The hit 
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count distribution is normalised so that the mean hit count, {N ohs } is unity. Uncorrelated Gaussian noise was then added to 
the T Q and U maps with a Q = a u = a T /V2 = 2.725/^ (iV obs ) (jK. 

We concentrate on the analysis of the E and B power spectra in this Section. With these parameters, the E and B mode 
power spectra become noise dominated at at I ^ 750 and i £ 50 respectively (see Figure 10). We apply an apodised KpO 
mask, as shown in Figure 4 but recomputed for the smaller pixels appropriate to these maps. In addition, a further weight 
factor was applied 

- ^ (57) 



as described in E04, where cx^ is the noise level for a pixel in the Q or U map with the mean hit count, and erf is the actual 
noise level in pixel i. The parameter e/ controls the weighting scheme so that e/ = corresponds to inverse variance weighting 
and e/ — » co corresponds to equal weight per pixel (apart from minor variations in the pixel areas £li associated with the 
pixelisation scheme). 

The results plotted in Figure 10 show equation (55) applied to two sets of PCL estimates, one with equal weight per pixel 
and one with £/ = 0.2. (The cross-covariances between these estimates can be computed via a straigtforward generalisation of 
equations (45) - (47).) The solid lines in Figure 10 show the theoretical dispersion computed from the diagonal components 
of the covariance matrix (56). As can be seen from the scatter of the points, these estimates are entirely reasonable. In fact, 
in this example, the variances of the PCL estimates at high multipoles are relatively insensitive to the weighting scheme. At 
high multipoles, the e/ =0.2 weighting produces slightly (~ 10%-20%) smaller errors compared to equal weight per pixel. 
However, in the signal dominated regime at low multipoles, the £/ = 0.2 weighting gives much larger errors. This behaviour 
is similar to the behaviour for the temperature anisotropics illustrated in Figure 14 of E04. 

To apply the QML estimator, the maps at high resolution were smoothed with a Gaussian beam of width 6 3 — 5° and 
repixelised onto a low resolution map with pixel size 9 C = 3.5°. The QML estimator (48) requires the noise covariance matrix 
Nij of the low resolution data. Since the noise in each high resolution map is diagonal, crfSij, the noise covariance matrix of 
a low resolution map is given by 

N H = XX (2 ^ + 1)(2£ ' + l)^(cose ip )P r (cose 3P )cxp (-^ 2 0i) cxp (-^' 2 #i) An '. ( 58 ) 

where p denotes the pixel number in the high resolution map, Q p is the solid angle of pixel p and 9^ is the Gaussian smoothing 
scale of the low resolution maps. In the small angle approximation, equation (58) simplifies to the easily computable expression 

*"5^E*«»(-£)-(-|j). < 59 » 

which is sufficiently accurate for the example discussed here. 

Figures 11a and lib show the PCL estimates at low multipoles from Figure 10 on a greatly expanded scale. The error 
bars on the points show the dispersions computed from the diagonal components of the covariance matrix The solid (red) lines 
show the E and B-mode power spectra for this particular simulation and the (green) dashed lines show the theoretical power 
spectra of the concordance ACDM model. The (purple) dashed line in Figure lib shows the theoretical power spectrum of 
the noise. 

Figures 11c and lid show analogous plots for the QML estimator, corrected for the low resolution smoothing 9l- In 
agreement with the tests shown in Figures 6 the errors on the QML estimates at low multipoles are considerably smaller 
than those for the PCL estimates. The covariance matrices for the QML estimates are very close to diagonal (c./. Figure 
7). Furthermore, the agreement between the E and B mode QML power spectrum estimates and the power spectra for this 
particular realisation (shown by the solid red lines), confirms the arguments given in Section 3.3 that the QML estimator 
will return almost the exact input power spectra even in the presence of a substantial sky cut. This is similar to the case of 
temperature power spectrum estimation, where a KpO-like sky cut can introduce a substantial 'estimator induced' variance 
(Efstathiou 2004b) at low multipoles comparable to the cosmic variance (c/. Table 1 of Efstathiou 2004b). For polarization, 
however, the problem of estimator induced variance is particularly acute because of E and B-mode mixing. If the B-mode 
amplitude is low, then the estimator induced variance of PCL estimators can easily overwhelm cosmic variance, whereas it 
is essentially negligible for a QML estimator. Notice that at multipoles £ ~ 15- 25 the PCL and QML estimates are closely 
similar (as expected from Section 3.2). At multipoles ^ ^ 25 the errors on the QML estimates become larger than those of the 
PCL estimator because of the large smoothing applied to the low resolution maps. However, there is quite a wide range of 
multipoles over which the estimates and errors are almost identical. The pixel size and smoothing scale of the low resolution 
simulations used in this example was chosen for computational speed so that the QML code could be run quickly on a single 
processor workstation. It would be entirely feasible to extend the QML calculation to multipoles of a few hundred or more 
using modern multi-processor computers. For many realistic experimental configurations this would provide a large range of 
overlap with PCL estimates at high multipoles. Thus, in practice, it would be extremely accurate simply to combine PCL 
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Figure 11. Figures (a) and (b) show the PCL estimates from Figure 10 at low multipoles on a greatly expanded scale. The error bars 
are computed from the diagonal components of the covariance matrix (56). Figures (c) and (d) show the QML estimates corrected for 
the additional smoothing applied to the low resolution maps from which they were computed. The solid (red) lines show the input E 
and B power spectra for this particular simulation, while the dashed (green) lines shows the power spectra of the fiducial ACDM model. 
The dashed (purple) line in Figure 11(b) shows the noise contribution to the B-mode spectrum. The noise contribution to the E mode 
spectrum is negligible at these multipoles. 



estimates using equation (55) and to overwrite the estimates at low multipoles, and associated components of the covariance 
matrices (56), with the corresponding QML estimates and covariance matrices. 

Alternatively, one can calculate the cross-covariances between the QML estimates and PCL estimates, and include the 
QML estimates as part of the data vector C r al used to compute the hybrid solution C( h . One then requires estimates of the 
cross-correlations between the PCL and QML estimates. These are quite complicated, but for completeness, we give analytic 
expressions here for the case of noise-free data. If we write the QML estimates in the form 



E 

Ve = 



QU 



UQ UU 
and define the quantities, 

U X = - ^(^(Mflft +iC QU (0 ip )R7')w i SU 



(60) 



(6l£ 
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= A Y,(C U Vip)KL + iC QV \B iv )R:*)wSli, (61b) 



= -4 ^(c Q («W™ ™ ^(M^Ha, (6ic) 



t# m = J2(c u (e ip )R^ - *c qc/ (M<>a ; (eid) 

then the cross-covariances are given by 

/Ar £ A f /\ - I L^m ZU plm U qlm l^ m \ U plm V qlm + U pl m V qlm) \ prt' ffiOn"! 



/A(^- B Ai, r \— ' I l^m ZV pim V qem 2~, m y V plm, U qem + V pim U qtm) prf' ffiOM 

As explained above, for most purposes, it should be sufficiently accurate simply to overwrite PCL estimates and covariance 
matrices at low multipoles with those computed from QML, avoiding the complexity of evaluating equations (62a) and (62b). 



6 CONCLUSIONS 

This paper generalises the analysis of power spectrum estimators presented in E04 to the case of polarisation. Simple analytic 
expressions involving the power spectrum of the square of the window function (the scalar approximation) are given for the 
covariance matrices of PCL estimates, including the effects of uncorrelated instrument noise. For noise free data, the scalar 
approximation is shown to give accurate error estimates for the _E-mode power spectrum at high multipoles for realistic 
(KpO-like) sky-cuts, and also for the B-mode power spectrum provided the Q and U maps are apodised appropriately and 
the amplitude of the B-mode is high enough. For noise dominated data, the scalar approximation povides extremely accurate 
covariance estimates for both the E and B modes at high multipoles. 

The results presented in Section 3 establish certain relationships between QML and PCL estimators, in particular, the 
statistical equivalence of QML and PCL estimators in the noise-free and noise-dominated limits. This analysis parallels the 
discussion of temperature power spectrum estimates given in E04. However, in the noise-dominated limit, an optimal PCL 
polarization estimator can only be constructed for the special case of uncorrelated noise and identical noise levels in the Q and 
U maps, (cr Q )j = (<J U )i • Fortunately, for a Planck-type experiment, both of these assumptions should be reasonably accurate 
at high multipoles. 

An analytic discussion of the noise properties of CMB temperature maps for a PZoncfc-like scanning strategy with realistic 
'1//' noise is given by Efstathiou (2005). For a Planck-type experiment, the noise should be accurately white at high multipoles 
(see also Stompor and White, 2004). Low frequency '1//' noise introduces striping in the maps with a characterstic spectrum 
that varies approximately as ACe oc l/£ at low multipoles. This low frequency behaviour for the temperature power spectra 
has been verified in many simulations of Planck- like experiments (e.g. Burigana et al. 1997, Maino et al. 1999; Keihanen et 
al. 2004; Efstathiou 2005). It is also seen in the X, E and B power spectra from full-scale simulations of the polarised Planck 
217GHz detectors (Ashdown et al. 2006). 

For the knee- frequencies 1 1 and noise levels expected of Planck detectors, the noise power spectra for both temperature 
and polarization will be accurately white-noise above multipoles of I Si 20. For the temperature power spectrum, the effect 
of residual striping errors on the power spectrum at multipoles £ 20 should be negligible, and so it should be an excellent 
approximation to treat the noise as pure white-noise. However, the simulations of Ashdown et al. (2005) show that the effects 
of striping errors on the E and B-mode power spectra at low multipoles will be comparable to (or will dominate) the intrinsic 
signal. These errors should therefore be folded into estimates of the power spectrum covariance matrices. This can be done 
by evaluating the QML estimator on low resolution maps, provided one can adequately approximate the noise covariance 
matrices Nij for these low resolution maps. In Section 5 we showed that it is easy to compute Nij for a low resolution map 
if the noise is strictly uncorrelated. However, the analogous problem for Planck-Yike correlated noise has not yet been solved 
and is currently under investigation. 

Since the noise for a Planck-like experiment is accurately white at some characteristic multipole, a hybrid estimator based 
on PCL estimates at high multipoles and QML estimates at low multipoles, is well motivated. Furthermore, by applying a 
QML estimator at low multipoles, it is possible to eliminate almost completely the troublesome effects of E and B mode mixing 
associated with sky cuts. We have previously argued (E04, Efstathiou 2004b) that maximum likelihood estimators should be 

II The characteristic frequency above which the detector noise is approximately white. 
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used to determine the temperature power spectrum at low multipoles on a cut sky to eliminate the large 'estimator-induced' 
variance inherent to PCL estimators. The problem of estimator-induced variance is even more acute for PCL estimators of 
the polarization power spectra (c/CC05). In particular, for the PCL estimators discussed here, mode-mixing places significant 
limits on the amplitude of the tensor-to-scalar ratio, r, that can be probed even in a noise-free experiment. For example, the 
tests described in Section 3.4 show that it is difficult to probe below r ~ 0.1 using a PCL estimator applied to maps with a 
KpO sky-cut. As Smith (2005) has shown, it is possible to define fast estimators based on unambiguous E and B modes on a 
cut sky. However, applying a QML estimator to low (or intermediate) resolution maps is also fast and is very close to optimal. 
In fact, if one wants to incorporate QML power spectrum estimation into a Monte-Carlo chain, then as equation (29) shows, 
the time taken to evaluate a QML power spectrum is no greater than the time taken to evaluate a fast spherical transform. 

Acknowledgements: I thank members of the Cambridge Planck Analysis Centre, especially Mark Ashdown and Anthony 
Challinor, for helpful discussions. 
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